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ABSTRACT 


Alternate  integral  forms  for  the  cumulative  probability  distribution  in  terms  of 
the  characteristic  function  are  given.  In  particular,  forms  that  can  utilize  a  fast 
Fourier  transform  (FFT)  algorithm  and  special  forms  for  one-sided  probability  den¬ 
sity  functions  are  derived,  For  a  special  class  of  discrete  random  variables,  all 
integral  evaluations  are  over  a  finite  range.  Some  computational  aspects  of  uti¬ 
lizing  the  FFT  are  discussed. 
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ALTERNATE  FORMS  AND  COMPUTATIONAL  CONSIDERATIONS  FOR 
NUMERICAL  EVALUATION  OF  CUMULATIVE  PROBABILITY  DISTRIBUTIONS 
DIRECTLY  FROM  CHARACTERISTIC  FUNCTIONS 


1.  INTRODUCTION 

A  recent  report  [l]  on  numerical  evaluation  of  cumulative  probability  dis¬ 
tribution  functions  directly  from  characteristic  functions  (QF)  gave  the  cumula¬ 
tive  distribution  functions  (CDF)  in  terms  of  a  single  integral  on  the  CF  for  both 
continuous  and  discrete  random  variables  (RV).  In  this  report  some  alternate 
forms  for  the  CDF  in  terms  of  the  CF  will  be  presented,  with  an  aim  toward 
more  accurate,  efficient,  and  expeditious  calculations.  For  the  motivation  of 
this  study  and  utility  of  the  results,  as  well  as  numerical  examples,  see  Ref¬ 
erence  1. 


2.  ANALYSIS 

This  section  is  composed  of  five  subsections.  In  the  first,  general  dis¬ 
tributions  are  considered;  in  the  second,  specialization  to  a  nonnegative  random 
variable  is  made.  In  both  subsections,  forms  that  utilize  a  fast  Fourier  trans¬ 
form  (FFT)  are  derived  and  their  applicability  is  discussed.  In  the  third  and 
fourth  subsections,  discrete  random  variables  are  considered.  The  former 
subsection  shows  that  the  distribution  function  can  be  evaluated  entirely  in 
terms  of  finite  integrals;  the  latter  subsection  specializes  to  nonnegative  dis¬ 
crete  random  variables.  The  fifth  subsection  treats  some  computational 
aspects  of  the  FFT. 


2.1  GENERAL  DISTRIBUTIONS 

Let  RV  x  have  probability  density  funotion  (PDF)  p(x)  and  CF  f({): 

f(f)  -  J dx  exp(i{x)  p(x)  ,  (1) 

P<x)  "  /**  exP(-i*x)  f<f)  ’  (2) 
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(Integrals  without  limits  are  over  the  real  axis  from  -«  to  «  .)  The 
CDF  Pr(X)  is  defined  as  the  probability  that  RV  x  is  less  than  or  equal  to  X. 
The  modified  distribution  function  (MDF)  P(X)  is  defined  equal  to  Pr(X)  at 
points  of  continuity,  but  it  takes  a  value  midway  between  limit  values  on  either 
side  of  a  discontinuity. 

The  MDF  P(X)  can  be  obtained  from  the  CF  by[lt  Eq.  (7),  or  2,  Eq. 
(4.14)] 


00 

1  fd£ 


Im  |f($)  exp(-i{X)  J, 


all  X  . 


(3) 


If  we  attempt  to  remove  the  Imaginary  operation  from  under  the  integral  sign, 
we  obtain  an  infinite  integral  since  f(0)  =  1.  However,  if  we  express 

f<*)  -  [f({)  -  a({)]  +  a(()  ,  (4) 

where  a(0)  =  1,  and  split  (3)  into  two  integrals,  we  oan  move  the  Imaginary 
operation  out  of  the  first  integral  in  (3).  One  particularly  useful  choice  for 
a((),  which  results  in  a  dosed  form  expression  for  the  second  integral  in  (3),  is 

a1(f)  »  exp[hi{  -  J  9Z  f2],  i  >  0  ,  (5) 

where4*  n  and  <rz  are  the  mean  and  variance  of  RV  x.  The  mean  and  vari¬ 
ance  are  available  from  f'(0)  and  f'(0),  if  these  quantities  can  be  evaluated;  if 
not,  the  method  to  be  described  is  still  applicable  with  arbitrary  constants  used 
for  n  and  a2.  When  (4)  and  (5)  are  substituted  into  (3),  there  results  [3,  Eq. 
3.896  4;  integrate  both  sides  with  respect  to  b] 


P (X)  -  ♦! 


m  - 


}• 


exp(-lfX)),  all  X 


(«) 


’Actually,  *  and  a2  could  bo  aaainod  arbitrary  values  in  tbo  torn  (Sit  thia 
particular  choice  gives  a  a#cond*ordar  fit  to  f(f)  at  tbo  origin. 
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where 


♦(y)s 


exp(-t2/2) 


(7) 


is  the  Gaussian  CDF. 

Equation  (6)  is  now  in  a  form  where  an  FFT  can  be  utilized  in  the  integra¬ 
tion  on  {  $ce  Subsection  2.5).  This  equation  is  exact;  we  are  not  making  a 
Gaussian  approximation  in  (6).  There  is  no  problem  in  the  integration  at  I  »  0 
because,  for  the  ohoice  uT  zr.i  „2  as  the  mean  and  variance  of  RVx, 


*<l)  -  MS)  2 
- 1 - -  0(0  as  t  —  0+  . 


(8) 


Also,  since  |  *i(S)[  decays  as  expf-v^fS/g),  the  decay  of  the  left  side  of  (8) 
far  large  {  will  often  depend  on  the  deoay  of  f(()A ;  this  decay  will  dictate  how 
far  the  integral  in  *6)  must  be  carried  out  for  specified  aoouracy  in  P(X). 

Other  ohoices  for  a(f )  are  possible  and  sometimes  recommended.  For 
example,  if  the  mean  and  variance  of  RV  x  do  not  exist  (e.g. ,  p(x)  ■  r”1 . 

(1  +  x2)~l,  all  x),  we  might  choose 


Sjd)  -  exp(-bf),  S  >  0  . 


(») 


To  best  match  f(f)  near  the  origin,  we  could  choose 

b  ■  -T(0+)  ■  | f ((H) |  (assuming  f(0+)  real)  .  (10) 

Then  by  substituting  (4)  and  (0)  into  0)  [3,  Eq.  3.041 1],  we  get 


P(K)  *  j  ♦  J  arotanpt/b)  -  J  to  j 


H I)  -  vo 


expf-KX), 


)• 


aUX. 

00 


For  the  choice  of  b  in  (10), 


f<{)  -  a2(|) 
- 


as  ( 


so  no  problem  in  integration  arises  at  the  origin.  We  must  be  able  to  evaluate 
f  (0+)  in  this  case  so  that  b  is  known,  and  it  must  be  real.  In  cases  where 
f 1  (04)  is  lot  known  or  Is  infinite  (See  Appendix  A,  for  example),  the  above 
methods  are  inapplicable,  and  special  techniques  such  as  subtracting  out  the 
singularity  are  required. 


2.2  NONNEGATIVE  DISTRIBUTIONS 

When  RV  x  is  limited  to  nomy»gative  values,  some  simplifications  in  the 
general  form  (3)  occur.  (The  case  of  nonpositive  RV  x  can  be  treated  in  a 
similar  fashion.)  First,  if  X  <  0  in  (3),  then  P(X) »  0.  Letting  X  »  -a  yields 

2  *  »  Cy**  •  >  0  »  f13) 


where  subscripts  r  and  i  denote  real  and  imaginary  parts,  respectively. 
Employing  (13)  in  (3)  for  X  >  0,  we  get 


!  />• 


<«  sin(fX),  X  >  0  , 


(f)  oosffX),  X  >  0 


Thus,  the  MDF  P(X)  can  be  evaluated  te cm  knowiedge  of  either  the  real  part 
or  tbs  Imaginary  part  of  the  CFfft).  For  X  *  0,  neither  (14)  nor  (15)  is 
necessarily  valid,  and  are  must  resort  to  (I). 


4 


There  are  computational  reasons  for  choosing  (14)  over  (15),  or  vice 
versa.  The  first  has  to  do  with  ease  of  calculating  fr(()  versus  1^(0 .  For 
example,  in  Appendix  A,  for  p(x)  =  2/V  (1  +  x2)”1  for  x  >  0,  we  find  that 
fr({)  is  a  simple  exponential,  whereas  fj(()  is  a  sum  of  exponential  integrals. 
Converse  examples,  where  fj(()  is  simpler  to  compute,  can  also  be  found. 

i 

The  second  reason  has  to  do  with  the  rate  of  decay  of  fr(£)  versus  fj(() . 

We  have 

fr($)  PW  oos(fx)  P eM  cos(fx)  ■  J  dx  pe(x)  exp(i(x)  ,  (16) 


m  p(x)  sin(fx)  ~J  dx  Pq(x)  sin((x)  -  i“*  j  dx  p q(x)  exp(i|x)  , 

0  (17) 

where  subscripts  e  and  o  denote  even  and  odd  parts,  respectively.  Now,  if 
p(0+)  >  0,  then  p0(x)  is  discontinuous  at  the  origin,  and  fj(()  decays  only  as 
I"1  far  large  f  .  An  example  is 

p<x)  -  e"*,  x  >  0;  f(|)  -  (1  -  iff1 . 


«r(« 


ae> 


In  (18),  fr({)  deosys  as  (■*  for  large  f ,  giving  rise  to  an  integral  in  (14) 
that  oan  be  terminated  earlier  than  the  one  in  (16).  On  the  other  hand,  consider 
that  p(0+)  *  0  and  that  ppt)  and  its  derivative  are  continuous  excopt  at  the 
origin,  but  p'(0+)  >  0.  Than  Po(x)  and  its  derivative  are  oonttnuous,  whereat 
Pgpt)  is  dleoontinuous.  lathis  ossa,  fr(()  decays  only  as  (**  for  large  f. 
An  aacampla  is 

p(x)  -  se"*,  x  >  0;  m  -  <1  -  t()“*  , 
yt>  •  (l  “  ^ )  (l  ^  **)“*»  fitt)  “  ♦  <*)'*  .  &»> 
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H'jre  (15)  could  be  terminated  earlier  than  (14). 

The  third  reason  has  to  do  with  the  region  of  X  of  interest.  For  large  X, 
where  P (X)  is  near  unity,  Eq.  (15),  in  the  form 


*  P(X) 


■if- 


■y  fjlf)  cos({X), 


X  >  0 


(20) 


is  to  be  recommended,  since  it  is  an  alternating  sum  of  small  quantities  and  re¬ 
tains  significance.  Equation  (34),  for  large  X,  is  an  alternating  sum  of  large 
quantities  and  I^scs  significance.  But  for  small  X,  Eq.  (14)  would  be  recom¬ 
mended. 

Equation  (15)  can  be  immediately  manipulated  into  a  form  where  an  FFT 
can  be  utilized.  Namely, 


P(X)  =  1  -  7  Re 


{  J  d|  (  exPH£x)> 


exp(-i*X)> ,  X  >  0  . 


(21) 


From  Appendix  A,  we  have  fi(()/(  ~  mx  as  (  — ■  0,  if  mx  exists  and  is  finite. 
If  we  attempt  to  express  (14)  in  the  form 


Im 


{»/? 


fr({)  exp(i£X) 


)• 


we  obtain  an  integral  that  does  not  converge  at  the  origin.  However,  if  we  ex¬ 
press 


yi>  58  [fr«)  -  *>«)]  +  b«)  , 


(22) 


where  b(0)  *#  1  and  b(£)  is 


real,  then  (14)  becomes 
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P(X) 


fun  j/< ctf 


fr(i)  -  b«> 


exp(i£X) 


>}•*/■ 


d*  b«) 


sin(£X> 

“ T“ ’ 


X  >  0 


(23) 


and  an  FFT  can  be  used  on  the  first  integral.  Preferably,  the  second  integral 
should  be  integrable  in  closed  form.  A  particularly  useful  choice  is 


b(£)  =  exp 


(  >  0  , 


(24) 


where*  m2  is  the  mean-square  value  of  RV  x.  This  quantity  is  available  from 
fJJ(O)  if  it  can  be  evaluated.  When  (24)  is  substituted  into  (23),  we  get  (See  (3) 
through  (7) ) 


P(X)  -  1  +  |ta  |/df  —  **P[  -  -  ^ 


l 


exp(i{X)>,  X  >  0  . 

(25) 


) 


The  function 


fr({)  -  exp 


tM. 


0  as  £  —  0+ 


in  (25)  if  the  mean-square  value  m2  exists,  hi  many  cases,  it  decays  as 
fr($)/£  for  large 

The  fact  that  MDF  P(X)  can  be  obtained  from  either  the  real  or  imaginary 
parts  of  the  CF  for  a  nonnegative  distribution  are  manifestations  of  the  fact  that 
fr(£)  and  fj(£)  can  be  found  from  each  other;  in  fact,  they  are  related  by  Hil¬ 
bert  transforms.  For  p(x)  =  0  for  x  <  0,  and  no  impulses  at  the  origin,  we 
see  that  [4,  p.  38] 


*A§  in  the  footnote  to  Eq.  (5),  M 2  could  be  assigned  any  convenient  value. 
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f({)  =  J 6x  p(x)  exp(i£x)  =  J  dx  p(x)  U(x)  exp(i$x)  =  3{p(x)  U(x)} 


=  3{pW)  •  3{U(x)(  =  m  •  [“  i«)  +  2^1  =  £  [f«)  +  » 

(26) 

where  U(x)  is  the  unit  step  function,  3  denotes  a  Fourier  transform,  « 
denotes  convolution,  and  H  denotes  a  Hilbert  transform.  Therefore, 

f«>  =  i*  {m\  , 

or 

=  *r«> -  -* )y*>l  •  w 

For  the  cases  when  p(x)  contains  an  impulse  at  the  origin  of  area  c0,  the 
first  part  of  (28)  is  still  correct,  but  the  second  part  is  incorrect  by  the  addi¬ 
tive  constant  cQ.  However,  we  can  still  find  fr(£)  from  fj({)  by  utilizing  the 
fact  that  fr(0)  =  1.  Thus,  either  the  real  or  the  imaginary  part  of  the  CF  con¬ 
stitutes  complete  knowledge  about  the  MDF  in  the  case  of  a  nonnegative  dis¬ 
tribution. 


2.3  DISCRETE  DISTRIBUTIONS 

In  this  subsection,  the  RV  x  is  restricted  to  take  on  values  that  are  mul¬ 
tiples  of  some  fundamental  increment  A,  and  can  be  either  positive  or  nega¬ 
tive.  Although  the  equations  in  Subsection  2 . 2  are  applicable  here,  it  is  advan¬ 
tageous  to  have  forms  for  the  distribution  function  that  require  finite  integrals 
rather  than  infinite  ones.  We  have  for  the  PDF 

p(x)  *  “  kA)  *  (29) 

k 
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where  the  sum  ranges  from  -*>  to  oo  .  The  CDF  Pr(M)  for  integer  M  is 
given  in  Reference  1,  Eqs.  (20)  through  (26).  All  the  integrals  are  finite  inte¬ 
grals  except  for  the  one  in  (2C)  for  MDF  value  P(0); 


P(0)  = 


We  now  rectify  this  situation  and  obtain  a  finite  integral  for  P(0)  also. 
From  Reference  1,  Eq.  (15), 


ck  =  2“  J  d*  f(t)  exp(-ikA{)  . 


(30) 


Therefore,  by  using  Appendix  B  and  f(-$)  *  f*(£),  we  get 


c 

o 


tan(Af/2)  ’ 


(31) 


Then,  we  obtain  the  desired  result 
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P<°>  ■  t  \+K  -l-Trf'* 


k  =  -< 


fj«) 

tan(A£/2) 


(32) 


As  £  —  0+,  the  integrand  of  (32)  approaches  2mx/A  if  /xx  exists  and  is 
finite.  (There  is  no  integral  expression  for  P(0)  in  terms  of  fr(£).  Since  fr(£) 
is  the  Fourier  transform  of  p  (x)  (See  Eq.  (16)),  and  since 

c 


/  dxpe«  =  i, 


irrespective  of  the  form  of  p(x),  fr(£)  contains  no  information  about  P(0). 
This  is  analogous  to  the  general  distribution  case  where  P(0)  follows  from 
(3)  as 


00 

»«  =  ;■;/ ae  yji/t .  > 

0 


2.4  NONNEGATIVE  DISCRETE  DISTRIBUTIONS 

When  RV  x  is  limited  to  nonnegative  values,  the  CDF  Pr(M)  takes  on 
forms  requiring  either  the  real  or  imaginary  parts  of  the  OF  for  its  evaluation, 
just  as  in  Subsection  2.2.  To  see  this,  we  note  that  ck  in  (23)  is  zero  for 
k  <  0.  By  letting  k=-m  in  (30),  we  get 

t/A  ff/A 

J  d£  sin(mA£)  f^f)  -  /  d£  cos(mA£)  ff(£)  for  m  >  0  .  (33) 

o  o 


When  we  employ  (33)  into  (30)  for  k  >  0 ,  we  get 


c 


k 


o 


k  >  0  , 


(34) 
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or 


ck“ 


2A 

jr 


/ 


r/A 


d|  sin(kA£)  f  ft),  k  >  0 


(35) 


Therefore,  the  CDF  Pr(M)  for  Integer  M  is  given  by 

t/L 

Pr(M)  =2  ck  =  #•/  d{  fr(i)[|  +  2  oos(kA{)] 
k=0  o  k=l 


where  we  have  used  (34),  (30),  the  fact  that  f(-£)  =  f*(f),  and  Eq.  1.342  2  in 
Reference  3.  Equation  (36)  enables  evaluating  the  CDF  in  terms  of  the  real 
part  of  the  CF  alone. 

To  represent  Pr(M)  in  terms  of  fj({),  \/e  first  note  that  for  nonnegative 
RV,  the  general  formula  for  P(0)  in  (32)  beoomes 


1 

2 


tan(Af/2) 


(37) 


Now, 


M 

Pr(M)  -  ^ok 
k-0 


“  o 

0 


<U  fjtt) 


sin(kA{) 
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(38) 


where  we  have  employed  (35),  (37),  and  Eq.  1.342  1  in  Reference  3.  Equation 
(38)  is  complementary  to  (36)  in  the  sense  that  only  the  imaginary  part  of  the 
CF  is  necessary  for  evaluating  the  CDF.  The  reasons  given  in  Subsection  2.2 
for  selecting  (36)  or  (38)  in  a  particular  application  are  again  relevant. 


2.5  USE  OF  FFT  FOR  FOURIER  TRANSFORMS 
Many  of  the  integrals  in  this  report  take  the  form 


/• 


g(t)  exp(-i2*ft)  . 


Suppose  a  limit  T  on  the  integration  can  be  found  such  that 


if. 


g(t)  exp(-i2rft) 


<  t  for  all  f  , 


(39) 


where  •  is  some  specified  tolerance  or  error.  Then,  attention  can  be  focused 
on  evaluating 


GT(f) 


'/* 


g(t)  exp(-i2irft) 


(40) 


Since  the  integration  in  (40)  is  over  an  interval  of  length  T,  it  is  seen  that 
Oj>(f)  will  undergo  a  significant  change  in  value  in  an  interval  no  smaller  than 
l/r  in  f.  Thus,  one  might  initially  anticipate  that  (40)  should  be  evaluated  at 
values  of  f  *  n/T,  n  •  1,  2, . . . .  However,  in  many  oases,  this  resolution, 
l/r i  may  be  much  too  fine,  and  be  the  result  of  satisfying  (39)  with  a  very 
small  « .  In  such  oases,  values  of  &r(f)  at  some  multiple  of  the  fundamental 
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resolution  may  be  satisfactory,  say  m/T,  where  m  is  an  integer.  Thus,  we 
might  be  interested  only  in  evaluating  G^n^p  )  ,  n  =  1,  2, -  But  from  (40) 


dt  g(t)  exp(-i2mmt/T) 


m-1  (k+l)T/m 

=  J  dt  g(t)  exp(-i2rnmt/T) .  (41 

k*0  kT/m 

In  mniring  the  substitution  u  =  t  -  kT/m  in  (41)  and  defining  the  collapsed 
function 


gcM  =  k=0 


]£  «(u  +k  m)’  02S  u  'T/mi 


,  otherwise 


we  note  that  (41)  becomes 


:(“  -  /  dtt  ®o(u)  “"('“'ffe  “)  • 


The  odlapsed  function  g_  is  obtained  from  g  by  "pre-aliasing"  g  into  the 
interval  T/m.  If  we  define  the  Fourier  transform  of  gQ  as 


G0» 


/m 

du  g  (u)  exp(-i2irfu)  , 


(44) 


then  (43)  yields 
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GT(nT)’Go("T)- 


In  terms  of  the 


(45) 


That  is,  G>p(f)  can  be  evaluated  at  f  3  m/T,  2m/T, ... 
Fourier  transform  of  the  collapsed  function. 


Now  suppose  g^u)  is  sampled  at  Increments  of  h  in  (44),  where 

rp  /_ 

h  =  ~  ,  and  weighting  |  wk|  applied  to  the  samples  in  an  effort  to  approxi¬ 
mate*  (44).  That  is, 

G0<f)  S  wk  go(kh)  exp(-12«fkh)  *  8o(f)  .  (46) 

k*0 

The  approximation  in  (46)  will  be  good  if  g^  and  the  exponential  are  sam¬ 
pled  frequently  enough.  Thus,  if  the  exponential  is  not  to  vary  by  more  than  a 
radian  between  samples,  we  require 


(47) 


When  (45)  through  (47)  are  combined,  the  desired  values  are  given  by 


h  ^  wk  go(kh)  exp(-i2rkn/M)  (48) 
k>0 


if 


< 


(48) 


By  defining 


*For  extnpla,  Stapaon'a  r«U  haa  a  1/3,  "Jb41  *  "2n*f8  *  *  1/3* 
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(50) 


*k 


hwk  gc(kh),  1  £  k  £  M  -  1 

b w  g (0)  +  bw  g  (T/m),  k  -  0 
0  0  MC 


we  can  express  (48)  as 


G,r(n  *  G0(n  y)  3  Gc(“  ?)  “  £  \  exP(“i2Tkn/M>  .  <51> 

which  is  an  M-point  discrete  Fourier  transform  (DFT)  of  the  sequence  jdjJ. 
The  factor  1/r  in  the  upper  bound  cm  |n|  in  (49)  is  due  to  the  aliasing  in  the 
frequency  domain  that  takes  place  at  |n|  ■  M/2.  In  fact,  letting  jw^}  be  the 
samples  of  waveform  w(t)  at  t  *  kh  and  W(f)  its  Fourier  transform,  it  can 
be  shown  that  (Appendix  C) 

G0(f)  «  W(f)  *y>0(<  -  kM  y)  •  (“) 

Thus,  the  value  fi0(y  ■y)  is  oomposed  of  at  least  two  overlapping  tails  of 
Gq(0.  In  order  to  avoid  this  aliasing,  we  must  observe  (49). 

To  summarise,  the  values  of  Gj>(f)  at  f  -  n$f  are  given  approximately 
as  an  M-point  DFT  in  0)1)  of  the  sequenoe  {dj.}  in  (50).  When  (42)  is  substi¬ 
tuted  into  (50),  this  sequenoe  can  be  expressed  as 


As  is  obvious  from  (59),  g(t)  must  still  be  evaluated  from  0  to  T  in  incre¬ 
ments  of  h,  that  is,  at  mM  ♦  1  values.  However,  collapsing  reduces  the 
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size  of  the  FFT  from  mM  to  M,  with  an  attendant  reduction  in  computation 
time  and  round-off  error.  This  method  is  related  to  one  given  in  Refer¬ 
ence  5,  |>.  81. 

In  applying  this  technique  to  numerical  integration  of  CF's,  since  the  ex¬ 
ponentials  take  the  form  exp(+ifX),  we  note  that  the  increment  in  X  at  which 
values  are  obtained  by  employing  an  FFT  are  2r/T,  or  2»m/T  for  coarser 
resolution  as  above. 


3.  CONCLUSIONS 

Several  alternate  forms  for  direct  numerical  evaluation  of  the  CDF  or 
MDF  from  the  CF  have  been  presented  that  have  utility  in  different  situations, 
including  ease  of  calculation,  rate  of  decay  of  the  integrands,  and  the  probabil¬ 
ity  region  of  interest.  Also,  the  speed  of  the  FFT  and  the  large  number  of 
values  of  the  distribution  functions  that  are  quickly  available  make  the  formulas 
presented  attractive  in  a  large  number  of  practical  applications. 

In  the  case  of  discrete  distributions  with  RV  that  can  take  on  positive  as 
well  as  negative  values,  all  integrals  for  the  CDF  are  finite  and  over  a  half¬ 
period  of  the  CF.  Reevaluation  of  the  sines  or  oosines,  as  in  (36)  or  (38)  for 
different  values  of  M,  can  be  avoided  if  one  notes  that 

sinJj^M  +  jJaJ  -  sinjj^M  -  ~^a  +  aj  «  sin  jj^M  -  “JaJ  oos[s]  +  cos  jj^M  -  |JaJ  sln[a], 

with  a  similar  rosult  for  cosine.  Thus,  if  s  table  of  sin  (a)  and  oos(a)  for  the 
values  of  a  (Af)  is  constructed,  this  recurrence  relation  can  be  used  to  obtain 
the  higher  order  M -dependence  required  in  (36)  and  (38)  without  reevaluating 
sines  and  cosines. 
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Then, 


*.(()  =  exp(-|{|)  , 

*,<(>  =  sgn«)  7  jexp(-|(|)  Ei(|(|)  -  exp(|(|)  Ei(-|*|)| 


Since 


Ei{~ If! )  -  In  1(1  +C  -  |(j  +“  |(|2  +0(|d3)  as  |(|  -  0, 
Ei(|d)  =  ln  |d  +  C  +  HI  +“|(|2  +0(|||3)  as  |d  “*  0  , 


there  follows 


sgn<()  III  (-In  |(|  +  1  -  C  +  “  |(|2)+  o(|||  3)j 
~7  *  as  |(|  -  0  . 


Therefore, 


fjtt)  cos((X)  2  /  j 


i-fe) 


as  Kl-0, 


which  is  unbounded,  but  into gr able.  So,  in  those  cases  where  nx  is  infinite, 
the  behavior  of  f[(()/(  at  tbs  origin  must  be  handled  carefully  in  order  to 
accurately  evaluate  the  integral.  One  possibility  is  to  subtract  out  the  singu¬ 
larity  and  integrate  it  analytically. 


^1.  S.  Gradshteya  and  I.  M.  Ryxbik.  Table  of  Integrate,  Series  and  Product*, 
Academic  Preaa.  New  York,  1965,  Eq.  (3.723  1  k 

A2Ibid. ,  Eqa.  (8.214  1)  and  (8.214  2). 
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Bi 

M.  J.  Lighthill,  An  Introduction  to  Fourier  Anslysis  rod  Generslited  Function* , 
Cambridge  University  Press,  New  York,  p,  21,  definition  t. 

B2Ibid,  p.  25. 

B^Ibid,  p.  60,  definition  22. 


^Ibid,  p.  66,  Theorem  26. 

®®Ibid ,  p.  48,  definition  19. 
®®Ibid,  p.  66,  Theorem  26,  Note. 
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-  In2,  n  =  0 


•,  n  *  0 


where  we  have  used  Eq.  (4.384  3)  by  Gradshteyn  and  Ryzhik.®7  Therefore,  the 
generalized  function  In  |  sin(x/2)  j  can  be  expressed  as 


In  sin  ®  = 


.IViE 

2  Lmd  n 


sgn(n)  inx 


n=-*> 
n  ^  0 


If  we  define  the  derivative  of  the  generalized  function  lnj  sin(x/2)|  as  the 
generalized  function  cot(x/2)/2,  differentiation  of  the  last  equation  yields  the 
expression  of  the  generalized  function  cot(x/2)/2  asB® 


2  00t(|)  ”  2  S  Sg“W  eil“  ' 


n=— oo 


(This  equation  says  that  the  spectrum  of  the  generalized  function  cot(x/2)  is 
the  odd  impulse  train. B9)  And  sinceB*® 


n=-oo 


we  obtain 


£  »<x  -  n2ir)  +  i  |  -  * +£ 


S.  Gradshteyn  and  I.  M.  Ryihik,  Table  of  Integral*,  Series  and  Products, 
Academic  Press,  New  York,  1965.  ” 

“Op.  cit.,  M,  J.  Li gh thill,  p.  28,  Theorem  15. 

®^Ibid.,  p.  66,  Theorem  26,  Eq.  (36). 

®*®Ibid.,  p.  67,  Example  38. 
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or 


oot(i) 

n=— oo 


in  the  sense  of  generalized  functions. 
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Appendix  C 


ALIASED  SPECTRUM 


If  we  define  the  infinite  impulse  train 

3h(t)  =^|a(t  -  nh) 
n 

and  use  the  time-limited  character  of  g^,  as  given  in  (42),  it  is  possible  to 
manipulate  (46)  as  follows: 

M 

fl0«  “  h  X]  wk  Sc(kh)  exp(-i2rfkh) 
k=0 

T 

= J  dt  w(t)  gc(t)  h  «h(t)  exp(-i2xft) 


*/dt  w(t)  gc(t)  h  ^(t)  exp(-i2rft) 


-  3  )w(t)  gQ(t)  h  ^(t)}  -  W(f)  •  Go(f)  •  il/h( f) 


■  »•  •£«,(* -ki)- 

Kc 

Using  h  ■  t  we  see  that  (52)  results. 
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